Spatial–Temporal Correlations between Soil pH and NPP of Grassland Ecosystems in the Yellow River Source Area, China

In recent years, ecological concerns such as vegetation destruction, permafrost deterioration, and river drying have been paid much more attention to on the Yellow River Basin in China. Soil pH is regarded to be the fundamental variable among soil properties for vegetation growth, while net primary productivity (NPP) is also an essential indicator to reflect the healthy growth of vegetation. Due to the limitation of on-site samples, the spatial–temporal variations in soil pH and NPP, as well as their intrinsic mechanisms, remain unknown, especially in the Yellow River source area, China. Therefore, it is imperative to investigate the coupling relationship between soil pH and NPP of the area. The study coupled MODIS reflectance data (MOD09A1) with on-site soil pH to estimate spatial–temporal variations in soil pH, explore the response of NPP to soil pH, and assess the extent to which they contribute to grassland ecosystems, thus helping to fill knowledge gaps. Results indicated that the surface spectral reflectance for seven bands could express the geographic pattern of soil pH by applying a multiple linear regression equation; NPP exhibited an increasing trend while soil pH was the contrary in summer from 2000 to 2021. In summer, NPP was negatively correlated with soil pH and there was a lag effect in the response of NPP to soil pH, revealing a correlation between temperate steppes > montane meadows > alpine meadows > swamps in different grassland ecosystems. In addition, contribution indices for temperate steppes and montane meadows were positive whereas they were negative for swamps and alpine meadows, which are apparent findings. The contribution index of montane and alpine meadows was greater than that of temperate steppes and swamps. The approach of the study can enable managers to easily identify and rehabilitate alkaline soil and provides an important reference and practical value for ecological restoration and sustainable development of grassland ecosystems in alpine regions.


Introduction
Vegetation plays an important role in the exchange of carbon, water, and energy between the soil and the atmosphere through surface albedo, roughness, and evapotranspiration [1]. Vegetation properties are influenced by a variety of environmental factors [2], especially soil properties related to climate regulation and adaptation, biodiversity conservation, water filtration, and carbon sequestration [3][4][5]. There are intricate interactions between vegetation and soil properties at different scales [6], which have become a hot topic for ecological research in the context of global change [7,8].

Remote Sensing Data Acquisition and Processing
The Moderate Resolution Imaging Spectroradiometer (MODIS) Terra surface reflectance products (MOD09A1) and MODIS NPP products (MOD17A2) were used in this study with a spatial resolution of 500 m and intervals of 8 days, which can be downloaded from the United States Geological Survey (USGS) (https://lpdaacsvc.cr.usgs.gov/appeears/explore) (accessed on 4 March 2022). The MOD09A1 product provided a systematic correction for atmospheric factors including aerosols, gases, and Rayleigh scattering to estimate the surface spectral reflectance of seven bands (referred as B1∼B7) [27,28].
The monthly data sets were constructed by adopting the maximum value composite method (MVC) for the original reflectance and NPP data from 2000 to 2021 to minimize the influence of atmospheric and cloud contamination [29]. The data from June to August were then averaged to synthesize the summer data. A 90 m DEM (SRTMGL3 DEM) was also downloaded from the USGS. The boundaries of the administrative area and the Yellow River basin were downloaded from the Resource and Environmental Science and Data Center (https://www.resdc.cn/) (accessed on 23 April 2022).

Soil Samples Data Acquisition and Processing
The soil type, grassland ecosystem type, landscape characteristics, etc., were fully considered in the sampling design. As the reflectance of soil properties and vegetation is more significant during the peak growing season, our sampling period was end of July to mid-August, and the distribution of sampling points is shown in Figure 1b. A total of 130 soil samples were collected from the 0-10 cm soil layer. Each sampling point was uniformly numbered, located using GPS, and latitude and longitude coordinates were recorded. Based on the resolution of MODIS data, two sample quadrats were selected to collect biomass and soil samples at a 500 m × 500 m sample point, and the average value of two sample quadrats spaced 100 m apart from each other was used as the observation

Remote Sensing Data Acquisition and Processing
The Moderate Resolution Imaging Spectroradiometer (MODIS) Terra surface reflectance products (MOD09A1) and MODIS NPP products (MOD17A2) were used in this study with a spatial resolution of 500 m and intervals of 8 days, which can be downloaded from the United States Geological Survey (USGS) (https://lpdaacsvc.cr.usgs.gov/ appeears/explore) (accessed on 4 March 2022). The MOD09A1 product provided a systematic correction for atmospheric factors including aerosols, gases, and Rayleigh scattering to estimate the surface spectral reflectance of seven bands (referred as B1∼B7) [27,28].
The monthly data sets were constructed by adopting the maximum value composite method (MVC) for the original reflectance and NPP data from 2000 to 2021 to minimize the influence of atmospheric and cloud contamination [29]. The data from June to August were then averaged to synthesize the summer data. A 90 m DEM (SRTMGL3 DEM) was also downloaded from the USGS. The boundaries of the administrative area and the Yellow River basin were downloaded from the Resource and Environmental Science and Data Center (https://www.resdc.cn/) (accessed on 23 April 2022).

Soil Samples Data Acquisition and Processing
The soil type, grassland ecosystem type, landscape characteristics, etc., were fully considered in the sampling design. As the reflectance of soil properties and vegetation is more significant during the peak growing season, our sampling period was end of July to mid-August, and the distribution of sampling points is shown in Figure 1b. A total of 130 soil samples were collected from the 0-10 cm soil layer. Each sampling point was uniformly numbered, located using GPS, and latitude and longitude coordinates were recorded. Based on the resolution of MODIS data, two sample quadrats were selected to collect biomass and soil samples at a 500 m × 500 m sample point, and the average value of two sample quadrats spaced 100 m apart from each other was used as the observation value for that sample point. Soil pH was estimated by preparing a suspension with a soil/water ratio of 1:5 [30] using an air-dried and 1 mm sieved soil sample and deionized water.

Regression Analysis
The Extract by Mask tool in ArcGIS 10.5 was used to generate the B1∼B7 and NPP raster images. Then, the Extract Values to Point tool was used to extract the point data of B1∼B7 and NPP from the raster images that correspond with soil pH point data. Outliers were removed from 130 samples and sorted by sampling points from highest to lowest, with 2/3 of the samples selected as the calibration set and 1/3 of the samples selected as the validation set at equal intervals for model building and accuracy verification, respectively. The normal distribution histogram showed the data suitable for regression analysis. To ensure the reliability of the model, the MODIS data used were kept close to the acquisition time of the field data.
Correlation analysis of the values of the B1∼B7 at the sampling points and the corresponding pH indicated that there was an internal relationship between the pH and the reflectance (Table 1). Curve estimation was carried out using B1∼B7 as the independent variable and pH as the dependent variable of the calibration set [31]. The primary model types for the soil pH are listed in Table 2, the performance of the models was compared and selected based on the coefficient of determination (R 2 ) and root mean square error (RMSE). The coefficient of determination (R 2 ) was used to evaluate the relationship between the observed and predicted values, while the root mean square error (RMSE) measured the degree of inaccuracy of prediction [32]. The formulas for R 2 and RMSE are as follows [33]:  In Equations (1) and (2), R 2 is the coefficient of determination, RMSE is the root mean square error, i is the number of the sample point, n is the total number of sample points, V oi is the observed value at sample point i, V pi is the predicted value at sample point i.
B1∼B7 data from the validation set were calculated using multiple linear regression model to derive the predicted soil pH, which was then regressed against the observed soil pH (Figure 2), demonstrating that prediction accuracy was achieved. B1∼B7 data from the validation set were calculated using multiple linear regression model to derive the predicted soil pH, which was then regressed against the observed soil pH (Figure 2), demonstrating that prediction accuracy was achieved.

Trend Analysis
In the analysis of the raster data from 2000 to 2021, the least squares method and the F-test based on MATLAB software can be used to obtain the trend and significance test for each grid by iterating image by image. A significant trend was observed when the regression coefficients passed the significance test (F-test, p < 0.05) [34].

Cross-Correlation Analysis
In this work, the spatiotemporal coupling between pH and NPP in the summer of 2000-2021 was represented using the Pearson correlation coefficient, and the time-lag correlation was applied to quantify the time-lag effect of pH on NPP [35,36]. Correlation and significance analysis of raster data on the time scale was performed using the MATLAB software package.

Contribution Index
Some researchers have used the contribution index in their research [37,38], which is employed in this paper to measure the degree of contribution of each grassland ecosystem type to pH and NPP. The formula for CI is shown in Equations (3) and (4).
In Equation (3), is the contribution index of grassland ecosystems to pH, is the average pH of the grassland type I, is the average pH of the entire study area, is the area of grassland type , is the area of the entire study area.
In Equation (4), is the contribution index of grassland ecosystems to NPP, is the average NPP of the grassland type i, is the average NPP of the entire study area, and as above.

Trend Analysis
In the analysis of the raster data from 2000 to 2021, the least squares method and the F-test based on MATLAB software can be used to obtain the trend and significance test for each grid by iterating image by image. A significant trend was observed when the regression coefficients passed the significance test (F-test, p < 0.05) [34].

Cross-Correlation Analysis
In this work, the spatiotemporal coupling between pH and NPP in the summer of 2000-2021 was represented using the Pearson correlation coefficient, and the time-lag correlation was applied to quantify the time-lag effect of pH on NPP [35,36]. Correlation and significance analysis of raster data on the time scale was performed using the MATLAB software package.

Contribution Index
Some researchers have used the contribution index in their research [37,38], which is employed in this paper to measure the degree of contribution of each grassland ecosystem type to pH and NPP. The formula for CI is shown in Equations (3) and (4).
In Equation (3), CI pH is the contribution index of grassland ecosystems to pH, pH i is the average pH of the grassland type I, pH Avg is the average pH of the entire study area, S i is the area of grassland type I, S is the area of the entire study area.
In Equation (4), CI NPP is the contribution index of grassland ecosystems to NPP, NPP i is the average NPP of the grassland type i, NPP Avg is the average NPP of the entire study area, S i and S as above.

Temporal Scale Changes in pH and NPP
In this work, we calculated and analyzed the temporal variations of pH and NPP at monthly, seasonal, and interannual scales using pixels as the minimum unit of calculation. As shown in Figure 3, NPP increases significantly in summer and from June to August (pink dashed line), while pH decreases more rapidly (yellow dashed line). The NPP showed periodic changes with significant fluctuations throughout the summer, with the maximum in 2020 and 2017, and the minimum in 2003 and 2009, respectively. In the different months of the summer, there were variations in the increase in NPP from June to August, with the rate in June and August being higher than that in July; the pH decrease in June and August was significantly higher than that in July, with similar fluctuating trends.

Temporal Scale Changes in pH and NPP
In this work, we calculated and analyzed the temporal variations of pH and NPP at monthly, seasonal, and interannual scales using pixels as the minimum unit of calculation. As shown in Figure 3, NPP increases significantly in summer and from June to August (pink dashed line), while pH decreases more rapidly (yellow dashed line). The NPP showed periodic changes with significant fluctuations throughout the summer, with the maximum in 2020 and 2017, and the minimum in 2003 and 2009, respectively. In the different months of the summer, there were variations in the increase in NPP from June to August, with the rate in June and August being higher than that in July; the pH decrease in June and August was significantly higher than that in July, with similar fluctuating trends.  Figure 4 shows the geographic patterns of mean pH and NPP in summer from 2000-2021. Throughout the study area, annual mean pH ( Figure 4a) roughly exhibited the distribution characteristics of a high value in the northeast and southwest, while the differences in other regions were on the contrary. The pH trend (Figure 4b) showed an increasing trend in a few northern parts, but a large portion of the territory showed a decreasing trend. The annual mean NPP was higher in the northeast (Figure 4c), indicating that the region has a large carbon sequestration capacity, while the opposite was true in the southwest. Furthermore, the NPP is increasing across the entire northeast region (Figure 4d).  Figure 4 shows the geographic patterns of mean pH and NPP in summer from 2000-2021. Throughout the study area, annual mean pH (Figure 4a) roughly exhibited the distribution characteristics of a high value in the northeast and southwest, while the differences in other regions were on the contrary. The pH trend (Figure 4b) showed an increasing trend in a few northern parts, but a large portion of the territory showed a decreasing trend. The annual mean NPP was higher in the northeast (Figure 4c), indicating that the region has a large carbon sequestration capacity, while the opposite was true in the southwest. Furthermore, the NPP is increasing across the entire northeast region (Figure 4d).

Correlation of pH and NPP
To measure the response of NPP to changes in soil pH, the correlation coef between them was calculated spatially based on pixels from 2000 to 2021 ( Figure 5 ing summer and from June to August, there were negative correlations in most re with some regions showing no correlation. The average correlation coefficient f whole study area was −0.407 (p < 0.05) during summer, indicating a moderate ne correlation (Table 3). Additionally, in terms of correlation by month, the strongest lation was found in June with a correlation coefficient of −0.408 (p < 0.05), follow August with a correlation coefficient of −0.331 (p < 0.05), and July with the lowest c tion of −0.317 (p < 0.05).

Correlation of pH and NPP
To measure the response of NPP to changes in soil pH, the correlation coefficient between them was calculated spatially based on pixels from 2000 to 2021 ( Figure 5). During summer and from June to August, there were negative correlations in most regions, with some regions showing no correlation. The average correlation coefficient for the whole study area was −0.407 (p < 0.05) during summer, indicating a moderate negative correlation (Table 3). Additionally, in terms of correlation by month, the strongest correlation was found in June with a correlation coefficient of −0.408 (p < 0.05), followed by August with a correlation coefficient of −0.331 (p < 0.05), and July with the lowest correlation of −0.317 (p < 0.05).   On the monthly scale, there was a negative time-lag correlation between pH and NPP ( Figure 6). From the point of view of positive time-lag correlation (Figure 6a-c and Table  4), all of which also showed a negative correlation. The correlation coefficient between pH in June and NPP in July was −0.462 (p < 0.05), followed by pH in June and NPP in August, which had a correlation coefficient of −0.242 (p < 0.05). The correlation coefficient between pH in July and August in July was weakest at −0.118 (p < 0.05). It can also be seen that pH in June and July had a greater effect on NPP in August, with a stronger effect in July.
The negative correlation coefficients with a time lag (Figure 6d-f and Table 4) show that there is a negative correlation between the current monthly NPP and the soil pH in the following month. The correlation coefficients between pH in July and NPP in June or pH in August and NPP in June were −0.317 and −0.203, respectively, and were weakly negative correlations, indicating that NPP in June was responsive to pH in July and August, i.e., the increase in NPP in June would lower pH in July and August. However, pH in August was a moderately negative correlation with NPP in July, which is probably because the study area belongs to the alpine region and reaches the most vigorous vegetation in July and August.  Table 4), all of which also showed a negative correlation. The correlation coefficient between pH in June and NPP in July was −0.462 (p < 0.05), followed by pH in June and NPP in August, which had a correlation coefficient of −0.242 (p < 0.05). The correlation coefficient between pH in July and August in July was weakest at −0.118 (p < 0.05). It can also be seen that pH in June and July had a greater effect on NPP in August, with a stronger effect in July.
The negative correlation coefficients with a time lag (Figure 6d-f and Table 4) show that there is a negative correlation between the current monthly NPP and the soil pH in the following month. The correlation coefficients between pH in July and NPP in June or pH in August and NPP in June were −0.317 and −0.203, respectively, and were weakly negative correlations, indicating that NPP in June was responsive to pH in July and August, i.e., the increase in NPP in June would lower pH in July and August. However, pH in August was a moderately negative correlation with NPP in July, which is probably because the study area belongs to the alpine region and reaches the most vigorous vegetation in July and August.

Temporal Variation of pH and NPP in Grassland Ecosystems
Trends in pH and NPP in the four grassland ecosystems have varied over the past 22 years (Figure 7). NPP increased in temperate steppes, montane meadows, and alpine meadows, while it slightly decreased in swamps (pink dashed line). The NPP of temperate steppes had the highest growth rate with a noticeable fluctuation regularity, while montane meadows and alpine meadows were a lower growth rate than temperate steppes. Furthermore, as shown in Figure 7, the trend of pH decreased in four grassland ecosystem types with swamps, montane meadows, and alpine meadows just about right, while temperate steppes were declining at a slower rate than the other types (yellow dotted line).

Temporal Variation of pH and NPP in Grassland Ecosystems
Trends in pH and NPP in the four grassland ecosystems have varied over the past 22 years (Figure 7). NPP increased in temperate steppes, montane meadows, and alpine meadows, while it slightly decreased in swamps (pink dashed line). The NPP of temperate steppes had the highest growth rate with a noticeable fluctuation regularity, while montane meadows and alpine meadows were a lower growth rate than temperate steppes. Furthermore, as shown in Figure 7, the trend of pH decreased in four grassland ecosystem types with swamps, montane meadows, and alpine meadows just about right, while temperate steppes were declining at a slower rate than the other types (yellow dotted line).

Correlation of pH and NPP in Grassland Ecosystems
To study the effects of pH on NPP in different grassland ecosystem types, the average correlation coefficients of different grassland ecosystems were extracted based on the aforementioned research results ( Figure 5). According to the results (Figure 8), the correlation coefficients between pH and NPP for the four grassland ecosystems exhibited negative correlations at seasonal scales. From the perspective of the significance (p < 0.05), the correlation of each grassland type was significant. In contrast to the moderate correlation

Correlation of pH and NPP in Grassland Ecosystems
To study the effects of pH on NPP in different grassland ecosystem types, the average correlation coefficients of different grassland ecosystems were extracted based on the aforementioned research results ( Figure 5). According to the results (Figure 8), the correlation coefficients between pH and NPP for the four grassland ecosystems exhibited negative correlations at seasonal scales. From the perspective of the significance (p < 0.05), the correlation of each grassland type was significant. In contrast to the moderate correlation for temperate steppes and montane meadows, followed by the weak correlation between pH and NPP of alpine meadows, and swamps with an extremely weak correlation in summer.

Contribution Index of Grassland Ecosystems to pH and NPP
By comparing the contribution indices of each grassland ecosystem type in summer (Figure 9), it could be seen that the contribution indices of pH and NPP were various in all grassland ecosystem types. As shown in (Figure 9), the contribution index of montane meadows and alpine meadows was higher than that of temperate steppes and swamps, as their proportions were 32.70 and 65.53, respectively. In addition, the obvious result was that the contribution indices of NPP for temperate steppes and montane meadows were positive whereas the swamps and alpine meadows were the opposite. In addition, the contribution indices of the pH for montane meadows and alpine meadows played an important role, in which mountain meadows was a negative effect on pH while that of alpine meadows was positive, implying that alpine meadows could lower pH and montane meadows were on the contrary. In terms of absolute value, the great influence of alpine meadows had the greatest impact on pH. Since the contribution index of pH for temperate steppes and swamps was extremely weak, it is not stated here.

Contribution Index of Grassland Ecosystems to pH and NPP
By comparing the contribution indices of each grassland ecosystem type in summer (Figure 9), it could be seen that the contribution indices of pH and NPP were various in all grassland ecosystem types. As shown in (Figure 9), the contribution index of montane meadows and alpine meadows was higher than that of temperate steppes and swamps, as their proportions were 32.70 and 65.53, respectively. In addition, the obvious result was that the contribution indices of NPP for temperate steppes and montane meadows were positive whereas the swamps and alpine meadows were the opposite. In addition, the contribution indices of the pH for montane meadows and alpine meadows played an important role, in which mountain meadows was a negative effect on pH while that of alpine meadows was positive, implying that alpine meadows could lower pH and montane meadows were on the contrary. In terms of absolute value, the great influence of alpine meadows had the greatest impact on pH. Since the contribution index of pH for temperate steppes and swamps was extremely weak, it is not stated here.

Modeling of pH and Its Spatial-Temporal Change
Rising temperatures caused by global warming have degraded the permafrost in the Yellow River source area, having a profound impact on the vegetation and soil conditions [39][40][41]. Because soil pH has a significant impact on plant growth [42], soil microbial changes [43], as well as soil nutrients [44], it has emerged as an important soil variable [23]. More importantly, recent studies have found that soil pH is an important indicator in grassland ecosystem studies [45]. Given the complexity of the ecological environment and the numerous challenges that researchers face when conducting field sampling, the method allowed spatial information to be derived from remote sensing analysis with a large amount of spatially continuous data [46,47].
In previous studies, salinity was estimated using remote sensing techniques [46,48], and these approaches were adopted in estimating pH in this study. The multiple linear regression equation is the superior model presented and it is feasible to estimate soil pH using it over the B1∼B7 (R 2 = 0.542 and RMSE = 0.473) ( Table 2). In other words, the model of soil pH was accurate based on the B1∼B7. The results were promising and showed that soil pH can be reasonably, and easily determined using the B1∼B7 images obtained from Terra surface reflectance products (MOD09A1).
The results of the multiple linear regression analysis between B1∼B7 and soil pH revealed a decreasing trend in pH during the growing season from 2000 to 2014 (Figure 3a). One possible reason was that the national project of returning grazing land to grassland had improved the vegetation cover in the study area since 2000 [49]. Additionally, the apparent decreasing trend in June and August reflected the strongest relationship between pH and NPP (Figure 3b,d), which was expected to lower the pH, implying that B1∼B7 was a better predictor of pH at the beginning and end of the growing season were better predictors of pH. The findings were similar to the research findings of [22] who reported that estimating soil pH at the beginning of the growing season was more effective in spring wheat fields. Furthermore, as indicated in our spatial visualization research (Figure 4a,b) the mean pH was high in the northeastern and southwestern parts of the

Modeling of pH and Its Spatial-Temporal Change
Rising temperatures caused by global warming have degraded the permafrost in the Yellow River source area, having a profound impact on the vegetation and soil conditions [39][40][41]. Because soil pH has a significant impact on plant growth [42], soil microbial changes [43], as well as soil nutrients [44], it has emerged as an important soil variable [23]. More importantly, recent studies have found that soil pH is an important indicator in grassland ecosystem studies [45]. Given the complexity of the ecological environment and the numerous challenges that researchers face when conducting field sampling, the method allowed spatial information to be derived from remote sensing analysis with a large amount of spatially continuous data [46,47].
In previous studies, salinity was estimated using remote sensing techniques [46,48], and these approaches were adopted in estimating pH in this study. The multiple linear regression equation is the superior model presented and it is feasible to estimate soil pH using it over the B1∼B7 (R 2 = 0.542 and RMSE = 0.473) ( Table 2). In other words, the model of soil pH was accurate based on the B1∼B7. The results were promising and showed that soil pH can be reasonably, and easily determined using the B1∼B7 images obtained from Terra surface reflectance products (MOD09A1).
The results of the multiple linear regression analysis between B1∼B7 and soil pH revealed a decreasing trend in pH during the growing season from 2000 to 2014 (Figure 3a). One possible reason was that the national project of returning grazing land to grassland had improved the vegetation cover in the study area since 2000 [49]. Additionally, the apparent decreasing trend in June and August reflected the strongest relationship between pH and NPP (Figure 3b,d), which was expected to lower the pH, implying that B1∼B7 was a better predictor of pH at the beginning and end of the growing season were better predictors of pH. The findings were similar to the research findings of [22] who reported that estimating soil pH at the beginning of the growing season was more effective in spring wheat fields. Furthermore, as indicated in our spatial visualization research (Figure 4a,b) the mean pH was high in the northeastern and southwestern parts of the study region between 2000 and 2021, while it increased the trend in a few northern parts, which may be related to geographic factors (Figure 1b), so further research was needed.

Impact of pH on NPP
With the development of 3S technology, the analysis of the driving forces of NPP in different regions and scales has become more extensive [50]. Temperature, precipitation [51], net radiation [52], and human activities [53,54] have all been shown to influence NPP. In particular, grazing as one of the human activities is a major factor for changes in NPP in grassland ecosystems [55]. Owing to the fact that the vegetation in the study area belongs to alpine grassland ecosystems, the temperature is the main limiting factor for vegetation growth [56,57]. Increasing temperature has alleviated temperature stress on vegetation and created favorable conditions for vegetation growth, which may have contributed to the significant increase in NPP over the past 20 years (Figures 3a and 4d).
In addition, the existing literature using both modeling and experiments has revealed that the effect of soil pH on NPP is quite significant [58], implying that soil properties are also major influencing factors on NPP. Our study of the correlation between soil pH and NPP addressed the problem on a regional scale with the help of remote sensing and filled a gap left by previous studies, which revealed a significant negative correlation between pH and NPP in summer from 2000 to 2021 (Figures 5 and 6). Other findings included that the correlation between pH and NPP was more pronounced in June and August than in July ( Figure 5, Table 3), due to the unique climatic characteristics of the study area, including the influence of high elevation and low temperature. In particular, June is the beginning of the growing season for grasslands and August is the end, while July belongs to the high growth season and may influence pH predictions (as described in Section 4.1), which may explain the result of the lower correlation between pH and NPP in July.
Similar to our results, soil pH had distinct effects on community stability in different grassland types [59] and the NPP of each grassland type varies substantially [60]. As illustrated in Figure 7, the decreasing pH in temperate steppes, montane meadows, and alpine meadows was accompanied by an increase in NPP over time, whereas swamps, on the contrary, decreased during 2000-2021. This indirectly implied that pH, as a soil property, was only one factor influencing changes in NPP. An increase in soil temperature [61], which can reduce soil moisture, increase evapotranspiration, and affect plant growth in swamps [62], was cited as the main reason for the change in grassland productivity, whereas temperate steppes, montane grasslands, and alpine grasslands were less susceptible to water stress. This could explain the simultaneous decline in pH and NPP in swamps. Moreover, the negative correlation was validated by the results of the correlation between pH and NPP ( Figure 8) and the contribution indices ( Figure 9) in various grassland ecosystems. The contribution index of montane and alpine meadows was higher than that of temperate steppes and swamps since the area of each grassland type played a leading role in calculating the contribution index ( Figure 9).

Limitations and Prospects
Although remote sensing analysis has become a standard method for monitoring vegetation and soils, the results were still ambiguous due to the complicated spatialtemporal variability of vegetation and soil biochemical processes [63], as well as the effects of resolution of remote sensing data [46]. This study used MODIS products with a moderate spatial resolution of 500 m to quantitatively assess the effects of soil pH on NPP at the time-series scale. The results will provide a scientific basis for the management, restoration, and sustainable development of grassland ecosystems in the Yellow River source area, China. However, as valuable as these methods are, they are still limited by the weak spatial resolution and the limited number of field survey samples.
Our results relate to broader spatial scales, and further studies should employ hyperspectral data for soil properties at smaller regional scales. Meanwhile, in future research, Gaussian process regression, random forest, artificial neural network [64], or kernel ridge regression [23] may be a better attempt to accurately estimate soil properties; data fusion methods should also be adopted to overcome these limitations [65], which would be a better guide for the study of the plant-soil relationship.
In addition, both NPP and pH vary seasonally [66][67][68][69] yet only summer was investigated, and there are various soil properties while just the pH response to NPP is studied in this paper. As a result, seasonal and monthly variations should be involved in future studies to explain more deeply the intrinsic mechanisms and drivers of coupled relationships in grassland ecosystems. Additionally, last but not least, the results of the study are only representative of our study area and further validation is required to verify that the results for the wider region are consistent with the current results.

Conclusions
The novelty of this study stemmed from the scarcity of previous research on the effects of soil properties on vegetation using remote sensing. In this work, the gap was bridged by studying the regional variability and spatial correlation of soil pH and NPP. The following results were obtained: (1) From 2000 to 2021, the annual mean pH decreased with increasing NPP in both time and space. (2) Soil pH and NPP showed a moderate negative correlation during the growing season, and their interaction revealed spatial heterogeneity.
(3) The correlation between pH and NPP was strongest in the temperate steppe, followed by montane meadows and alpine meadows, with swamps being the weakest in summer. The correlation of each grassland type was significant (p < 0.05). (4) The contribution indices of pH and NPP for each grassland ecosystem type were distinct, with the larger contribution indices of montane meadows and alpine meadows implying a more significant contribution for the entire study area. In addition, the evident outcome was that whilst the contribution indices for alpine meadows and swamps were negative, those for temperate steppes and montane meadows were positive.
The results demonstrated the advantages of remote sensing applications in soil and vegetation research, allowing rapid analysis of their spatial distribution and changes over time. In summary, the results of the technique can be expressed well in topsoil studies, and further research should be conducted to investigate the relationship between soil and vegetation in deeper soils. Furthermore, we should focus on grassland ecosystem types that are sensitive to soil properties and restore areas of poor vegetation cover where possible in future practical work. The cost-effective strategy allows land managers to save time and effort to restore vegetation health and should be widely used in future research.

Data Availability Statement:
The data that support the findings of this study are available from the corresponding author upon reasonable request.